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We discovered a new type of spiral wave solutions in reaction-diffusion systems — spike spiral 
wave, which significantly differs from spiral waves observed in FitzHugh-Nagumo-type models. We 
present an asymptotic theory of these waves in Gray-Scott model. We derive the kinematic relations 
describing the shape of this spiral and find the dependence of its main parameters on the control 
parameters. The theory does not rely on the specific features of Gray-Scott model and thus is 
expected to be applicable to a broad range of reaction-diffusion systems. 
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Formation of rotating spiral waves (rotors) is one of 
the most vivid and ubiquitous phenomena of nonlinear 
physics IMl- These waves are observed in nonlin¬ 
ear optical media |j^, chemical reactions of Belousov- 
Zhabotinsky type, catalytic reactions on crystal surfaces 
110 , and in a variety of biological systems: social amoe¬ 
bae Dictyostelium discoideum , Xenopus oocytes Q , 
chicken retina ||l5|] , and the heart of animals and man, 
where the formation of spiral waves is responsible for 
cardiac arrythmias and the life-threatening condition of 
ventricular fibrillation 00 - 

A generic model used to describe spiral waves is a pair 
of reaction-diffusion equations of activator-inhibitor type 

m 


Bf) 

T9-^ = - q{0,T],A), 

Tr,^ = L'^Ar] - Q{e,r],A), 


( 1 ) 

( 2 ) 


where 0 is the activator, i.e., the variable with respect 
to which there is a positive feedback; 77 is the inhibitor, 
i.e., the variable with respect to which there is a neg¬ 
ative feedback and which controls activator’s growth; q 
and Q are certain nonlinear functions representing the 
activation and the inhibition processes; I and L are the 
characteristic length scales, and rg and t,; are the char¬ 
acteristic time scales of the activator and the inhibitor, 
respectively; and A is the bifurcation parameter. A con¬ 
siderable amount of studies was done on the excitable sys¬ 
tems with FitzHugh-Nagumo-type kinetics (N-systems) 
(see, for example, and references therein). These 
systems are described by Eqs. (Q) and (||) with L = 0, 
and the nonlinearity in q such that the nullcline of Eq. 
(||) is N- or inverted N-shaped. Theory of the spiral waves 
in N-systems with a = re/r,, 1 was recently developed 
by Karma |p^ , p^ . 


The existence of spiral waves in excitable N-systems 
is due to the ability of such systems to sustain traveling 
waves, the simplest of which is a solitary wave — trav¬ 
eling autosoliton (AS) p|-p^. In N-systems the equation 
q{9,r],A) = 0 has three roots: On, 9^, and 9^, for fixed 
A and rj = rji. Eor a <§; 1 an AS consists of a front, which 
is a wave of switching from the stable homogeneous state 
9 = 9h and rj = iju to the state with 9 = ^max (Oh = On 
and 0 max = 0 i 3 for rji = rjh) whose width is of order I, 
and a back of width of order I that follows the front some 
distance w ^ I behind the front. Thus, in the AS the 
distribution of 0 is a broad pulse, while the value of rj 
slowly varies from r] = rjh to some value rj = rjm in the 
back of the pulse, and then slowly recovers from rjm to 
r]h behind it. In the limit a —> 0 the amplitude of the 
wave (the value of 0 max) becomes independent of a and 
the speed c does not exceed the value of order l/rg, with 
both 0 niax and c determined only by the nonlinearity in 

q- 

At the same time, many excitable systems are de¬ 
scribed by Eqs. (||) and (^ in which the nullcline of 
Eq. (0 is A- or V-shaped. Examples of such A-systems 
are the well-known Brusselator and Gray-Scott models 
of autocatalytic reactions and an example of a V-wstem 
is the Gierer-Meinhardt model of morphogenesis 0 . Re¬ 
cently, we showed for the excitable Brusselator [M and 
Gray-Scott model 0 that they are also capable of prop¬ 
agating traveling waves — traveling spike AS. The prop¬ 
erties of these AS are essentially different from those in 
N-systems. The distribution of 9 in an AS has the form of 
a narrow spike whose amplitude grows as a decreases and 
can become huge as a —> 0. In contrast to N-systems, in 
the spike r\ varies abruptly and then slowly recovers back 
to r\h far behind the spike. The speed of such an AS is 
always greater than Z/rg and goes to infinity as a —> 0 . 
Also, it is important to emphasize that in an AS the front 
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and the back are not separated by a large distance, as in 
N-systems. In this Letter we will show that in excitable 
A- or V-systems it may also be possible to excite steadily 
rotating spiral waves and develop a theory of such waves 
in the case a 1. 

To be specific, we will consider excitable {L = 0) Gray- 
Scott model, which is described by the following equa¬ 
tions 10 


^ - e, 

at 

+ l - 7 ?, 


(3) 
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where we chose I and rg as the units of time and space, 
respectively. Recently, we showed numerically that a 
steadily rotating spiral waves may be excited in this 
model at sufficiently small a From these simula¬ 

tions one can see that the spiral wave is a steadily ro¬ 
tating slightly curved sharp spike front of width of order 
1 which in the cross-section looks like a periodic travel¬ 
ing wave train (sufficiently far from the core). We would 
like to emphasize that it is this kind of the concentra¬ 
tion profiles that is experimentally observed in Belousov- 
Zhabotinsky reaction 

Let us derive the equation of motion for the traveling 
sharp front in Gray-Scott model with a ^ 1 . Since 77 
varies slowly outside the front, it can be replaced by a 
constant value t} = rji ahead of the front. Let us introduce 
self-similar variable z = p — ct, where p is the coordinate 
along the normal direction to the front. For definiteness 
we will assume that c > 0, what means that the front is 
moving in the +z direction. In the sharp front the value 
of 9 is large , so one can neglect the last two terms in 
Eq. (® . Therefore, in the presence of small curvature K 
Eqs. and ( 0 ) can be written as 


^ -I- (c -I- -I- AO'^'q - 6» = 0 , 

dz^ dz 

dz 


(5) 

( 6 ) 


These equations have to be supplemented by the bound¬ 
ary conditions 0 (±oo) = 0 and 77(-|-oo) = iji, where the 
infinity actually means sufficiently far ahead of the front 
compared to the front thickness. 

Let us introduce new variables 


e = a^/^e(^-±^Y\ f, = ^, 

\ c J rji 

and the following quantities 

(jh) 

In these variables Eqs. i) and (i become 


c = c + K. 


(7) 


( 8 ) 


dH ^d9 --2. „ 

+ c— + A 9 ^r] -9 = 0 , 
dz‘^ dz 

dz 


( 9 ) 

( 10 ) 


with the boundary condition ?7(-|-oo) = 1. Observe that 
all a-, A- and ?7i-dependences now enter only via the 
parameter A. 

Equations (^ and have been analyzed by us in 
0- There we showed that these equations indeed admit 
solutions in the form of a traveling spike. In the spike 77 
varies from 77 = 1 at 2; = -|-oo (ahead of the front) to some 
7 ? = 77niin at z = —oo (behind the front). The value of 
0 ^ 1 in the spike, so in the original variables we indeed 
have 9 ^ ^ 1 [see Eq. (j^)]. The speed of the 

front c as a function of A obtained from the numerical 
solution of Eqs. (^ and (10) is presented in Fig. [| 0 . 
One can see that this solution exists only for A > Ag, 
and its speed is always greater than c = Cmin, where 

Ag = 3.76, c„,in = 1.5 (11) 

The analysis of Eqs. (^) and also shows that 
= 0.05 at A = Ag and rapidly decreases 


Vn 


= Vu 


as A increases |19 . 

From Eq. (|8|) immediately follows that for small K the 
correction 5 c to the velocity c due to curvature is 


r ry' (1 A dc 

6 c = —K 1 H-^ 

I 25 dA 


( 12 ) 


where in the right-hand side A, c, and dc/dA are evalu¬ 
ated at AT = 0 . Also, as we showed in |j^, for A not in the 
immediate vicinity of Ag with good accuracy c = 0 . 86 A 
and 77 „iin = 0 . Then, going back to the original variables, 
we may write 


where 


c = Coo - aK, 


Coo = 0.86Aa a = 


(13) 


(14) 


Behind the sharp front the value of ndrops from 77^ to 


77 arc 

77 niin and 9 goes exponentially to zero . On the much 
longer time scale a~^ the value of 77 recovers from 77 i„in 
according to Eq. (||), in which outside of the front the 
term — 9 '^rj can be neglected. From this we obtain that 
after the front passed a point x at time U = ti{x) we have 


r]{x,t) = 1 - [1 - 77min(ii)]e 




(15) 


where 77niin = rjiflmin- In a steadily rotating spiral wave 
we must have r]{x,ti + T) = rii{ti) = const, where T = 
211 joj and ui is the angular frequency of the rotation of 
the spiral. Therefore, the spiral should be described by 
Eq. 0) with Coo = const, which is in turn related to ui. 
This equation was first analyzed by Burton, Cabrera, and 
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Frank (BCF) for growth of screw dislocations on crystal 
surfaces (see also §). They calculated the shape of 
the spiral and its frequency in the case when the tip of the 
spiral is at rest. Applying their results to Eq. we 

obtain oj = Qrjf, where for simplicity we used 
the expressions in Eq. (|l4[) and put r/^in = 0. Since 
Arji > in order for the front to be able to propagate, 
we must have w ^ 1, so one can expand the exponential 
in Eq. (^, and obtain iji = and w = 

The spatial step h of the spiral far from the 
core will be ft. = Notice that a similar 

method was recently used to analyze asymptotically the 
spiral waves in N-systems [ p’blJT^ . 

Comparing the results obtained above with Eq. (0, 
one can see that in order for the solution in the form of 
the traveling front to exist, one should have A > ^ 

1. On the other hand, for A ^ we have a; 3> 1, 

so ft will not have enough time to decay behind the wave 
front. This means that this kind of spiral wave may exist 
only at A ~ Notice that according to Eq. dO) 

we have Coo ~ 1 and ft ~ 1 for these values of A, so the 
formulas obtained above for the frequency of the spiral 
should only be correct qualitatively. 

These facts may seem rather surprising since they pre¬ 
dict that the spiral waves should exist only in a narrow 
range of the values of A far from the excitability thresh¬ 
old A = AbT = 3.76a^/^, which is obtained from Eq. 

for rji = rjh = ^ [0- Also, numerical simulations 
show that spiral waves exist in a wide range of the values 
of A. What we will show below is that the spiral waves 
actually exist for all values of A;, < A < . 

The thing is that in addition to the spiral wave solution 
whose tip is at rest, there may also be a solution whose tip 
moves along a circle of some radius i?, which can be large 
for Aa^/^ <C 1. The reason the front will not propagate 
inside the circle is that the tip is right at the propagation 
threshold. This means that for A we have r\i = 

ffl = 3.76a^/^A“^ in the limit a —> 0 (for R ^ 1 the 
corrections due to curvature can be neglected), so the 
frequency to is determined by Eq. ( P^ ) with rji = ri\. In 
particular, for ^ A we can expand the 

exponential and obtain asymptotically 

OJ = l.lQa^/'^A. (16) 

This equation shows that the value of uj lies in the range 
^ u; < 1, as should be expected. For large values of 
R the speed of the front far away from the core should 
only slightly exceed Cmin, so the step of the spiral will be 
ft = 5.4 q;“^/^A“^. Note, however, that because of the 
closeness to the threshold point A = A^ the expansion 
in Eq. is no longer justified and therefore the BCF 
theory, as well as Eq. (|l^), is not applicable to the spiral 
waves in this parameter range. This theory can be mod¬ 
ified by noting that close to A{, we have approximately 

/ A \ 

c = c„,in + ft A - Afe - , (17) 

\ ^^min / 


where 6 is a constant and the tilde quantities in the right- 
hand side are evaluated at AT = 0. The analysis of Eqs. 
(^ and ( p^ shows that ft = 1.4. Observe that this equa¬ 
tion reduces to the form of Eq. (H) only very far from 
the core. 

Following j2^, we rewrite Eq. ®) for the steadily ro¬ 
tating spiral in terms of the angle (j) between the tangent 
vector to the front and the radius vector as a function of 
the distance r to the origin. A straightforward calcula¬ 
tion shows that in these variables Eq. (0) becomes 


d(f> 

dr 


1 

— tan ( 
r 


+ - 


2c„ 


Afcft^ cos (j) - 


b‘^(A - Ab) - (Cniin - Ojr cos (jjf 


(18) 


with the boundary conditions (/)(-|-oo) = 7r/2 |23| and 
(piyR) = 0. The latter condition says that at its tip the 
front is normal to the circle along which it rotates, what 
follows from the physical considerations. Since the front 
at r = i? is at its propagation threshold, its normal ve¬ 
locity there should be equal to Cmin- This, together with 
the boundary condition at r = i? immediately gives us 
R = Cinin/w. Knowing the value of i?, one can then cal¬ 
culate A — Ab and find the value of uj for which it agrees 
with Eq. (p^. This will determine a unique value of uj. 
Numerical solution of Eq. (pq ) shows that for a; <g; 1 we 
have A —Aft ^ 1, so Eq. (|16[) should indeed be recovered 
in the limit a ^ 0 with <C A <C and the 

spiral wave solution is an involute of a circle of radius R. 

The solution of Eq. (hi) for uj = 0.29 is presented in 


Fig. I For this value of oj we found A — Af, = 0.34, 
which gives rji = 0.86, within 4% in agreement with Eq. 

Comparing these quantities with the results of the 
numerical simulations of Eqs. (^) and (^) for a = 0.1 and 
A = 1.5, in which this value of oj was observed, we find 
that the value of r]i agrees with the predicted one within 
3% accuracy. The speed Coo = 2.3 obtained from Eq. 
0) also agrees with that observed numerically within a 
few per cent accuracy. 

This agreement is quite remarkable considering the fact 
that at these parameters the spiral wave already under¬ 
went meandering instability. In fact, according to our 
analysis, steady rotation of the spiral requires fine-tuning 
of the value of rji at the tip of the spiral. Note that the tip 
of the spiral is not described by the interfacial equations 
derived above and thus is a rather singular object capa¬ 
ble of sudden movements on the smallest length scale. 
Thus, it is natural to expect that the tip trajectory in 
a meandering spiral may be rather abrupt. Notice that 
similar situation is observed in the simulations of models 
of cardiac tissue (see, for example, p^). 

In conclusion, we developed a theory of spike spiral 
waves in Gray-Scott model. Spike traveling waves are 
observed in a variety of excitable systems including nerve 
and cardiac tissue. Even though we performed an analy¬ 
sis of a concrete system, Eqs. 0 and dH) have general 
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character and thus are expected to apply to other A- and 
V-systems (see also ||9|,|l0|) and other excitable systems 
of different nature in which spike traveling waves are ob¬ 
served. Also, such waves can be expected in combustion 
systems and Belousov-Zhabotinsky reaction in continu¬ 
ous flow reactors. The thing is that although in Eqs. 
(H) and (||) describing these systems the activator null- 
cline may formally be N-shaped, for typical parameters 
the value of 0max may be several orders of magnitude 
greater than 6h , so the system effectively behaves as a A- 
or V-system. In particular, this is true for the models of 
systems with uniformly generated combustion material 
P,p^ and the two-parameter version of Oregonator . 
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A 

FIG. 1. Dependence c(A) obtained from the numerical so¬ 
lution of Eqs. and (|i^. 



X 

FIG. 2. Steadily rotating spiral wave. Results of the nu¬ 
merical solution of Eq. ( [l^ ) with u! = 0.29 and A — Ah = 0.34. 
The circle shows the core of the spiral. 
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